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ABSTRACT 


This  paper  presents  the  results  of  empirically  testing  8 alternative 

31 

multipliers  for  a multiplicative  congruentlal  generator  with  modulus  2 -1. 

c 

The  LLRANOOM  random  number  package^^W^  uses  one  of  the  multipliers,  the 
simulation  programming  language  SIMSCRIPT  II  uses  a second  and  the  remaining 
six  are  the  best  of  50  candidate  multipliers  studied  by  Hoaglln  (1976)  using 
the  theoretical  spectral  and  lattice  tests.  The  battery  of  tests  fall  to 
detect  any  departures  from  randomness  for  3 of  the  multipliers,  even  at  a 
0.20  significance  level.  This  group  Includes  the  multiplier  that  SIHSCRIPT  II 
employs.  However,  another  of  the  3 superior  performers,  397204094,  requires 
only  78  percent  of  the  computing  time  that  the  SIMSCRIPT  II  multiplier  does 
and  is  the  second  most  efficient  computationally  of  all  8 multipliers. 


1.  Introduction 
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This  paper  presents  the  results  of  empirically  testing  8 multipliers 
suggested  In  the  literature  for  use  In  a multiplicative  congruentlal  pseudo- 
random number  generator  with  modulus  Generators  of  this  type  are  In 

common  use.  although  Marsaglla  (1968)  has  shown  that  all  such  congruentlal 

Q 

generators,  whether  they  be  of  modulus  2 or  of  prime  modulus, 
possess  flaws  that  make  their  theoretical  properties  differ  from  those  of 
an  Ideal  source  of  random  numbers.  Since  these  departures  from  Ideal  properties 
conceivably  could  cause  serious  errors  In  practice,  submitting  these  generators 
to  empirical  testing  provides  a way  of  evaluating  their  performance. 

Consider  the  linear  congruentlal  multiplicative  generator 

(1 ) E AZ^  1 (mod  M) 

31 

where  M ■ 2 -1  . In  order  that  (1)  generate  all  Integers  In  [1.  M-1] 

before  cycling,  A must  be  a primitive  root  of  N [ 8 ].  Although 

31 

{V^  ■ Z^/M  ; 1 ■ 1,...,M  - 1}  denotes  a sequence  of  2 - 2 distinct 

-31 

fractions  In  (0,1).  with  spacing  2 and  each  with  a finite  binary 
representation,  this  sequence  Is  not  the  one  encountered  In  practice. 

Because  of  the  need  to  assign  a byte  to  the  exponent  of  a floating  point 
number,  a common  procedure  on  IBH  360/370  system  computers  generates  the 
fractions  (U^  ■ (2[Z^/2®]  + 1 )/2^^  ; i « l,...,M  - 1}  where  [e]  denotes  the 

Integer  part  of  e . Then  the  sequence  {U^}  has  2 -2  fractions  per  cycle 

^Sce  Learmonth  and  Lewis  (1969)  and  Payne,  et  a1.  (1969). 
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that  assume  2^^  values  in  [1/2^^,  1 - 1/2^^]  in  increments  of 

1/2^^  = 0.119209...  X 10’^.  Moreover,  the  fractions  1/2^^  and  1-1/2^^  occur  2®-l 

O 

times  per  cycle  whereas  each  of  the  other  fractions  occurs  2 times  per 
cycle,  indicating  only  a minute  departure  from  uniformity.  This  density 
of  points  seems  sufficient  for  most  purposes. 

Although  the  density  consideration  in  (0,1)  is  important,  the  issue  of 
randomness  is  paramount.  Presumably  one  would  like  to  choose  a multiplier  A 
such  that  treating  {U^}  as  a sequence  of  i.i.d.  random  variables  from  U(0,1 ) 
introduces  incidental  error.  Lewis,  et  a1.  (1969)  recommend  A = 16807  and 

fail  to  detect  departures  from  the  assumptions  of  independence  and  unifoi  ity 
in  their  empirical  testing.  Learmonth  and  Lewis  (1973)  use  this  multiplier 
in  their  random  number  generator  LLRANDOM,  as  do  the  discrete  event  simulation 
programming  language  SIMPL/1  [7]  and  APL  [9].  Payne,  et  a1.  (1969)  recommend 

A « 630360016  and  claim  that  their  testing  shows  no  departures  from  assumptions. 

Recently  Hoaglin  (1976)  screened  50  primitive  roots  of  M using  the 
spectral  [3]  and  lattice  [15]  tests.  These  theoretical  tests  provide  an 
indication  of  the  relative  desirability  of  alternative  multipliers  for  gen- 
erating k-tuples.  The  present  study  reports  on  how  the  Lewis,  et  al.  choice 
A ■ 16807  (multiplier  I),  the  Payne,  et  al.  choice  A - 630360016  (multiplier  II) 
and  the  best  6 (multipliers  III  through  VIII)  chosen  from  Hoaglin’s  study  fare 
when  subjected  to  identical  empirical  testing. 

2.  Testing  Procedure 

I 

Three  hypotheses  were  tested:  | 

• {U^}  Is  a sequence  of  i.i.d.  random  variables.  ^ 

H2  • has  a uniform  distribution  on  (0,1). 

! 

- 


Hj.  H has  a uniform  distribution  on  the  unit 

square. 

For  each  multiplier  the  data  base  consisted  of  n = 100  nonoverlapping 
samples  or  replications  each  of  N = 200,000  . For  each  of  the  100 

replications  and  each  hypothesis  a statistic,  whose  asymptotic  distribution 
was  known,  was  computed  from  the  200,000  observations.  Then  the  100 
statistics  for  a given  generator  and  specific  hypothesis  were  subjected  to 
a battery  of  goodness-of-fit  tests  designed  to  detect  departures  of  their 
empirical  cumulative  distribution  functions  (cdf's)  from  their  corresponding 
theoretical  cdf's. 

3.  Runs  Up  and  Down  Test  Statistic 

To  test  we  relied  on  runs  up  and  down  statistics.  Let  denote 
the  number  of  runs  up  and  down  of  length  1 for  1 ^ 1,...,6  . Then  under 
the  quanti ty 

(2)  R-  I c^.[R^  - E(R^)][Rj  - E(R  )] 

1,j=1 

asymptotically  has  the  chi-square  distribution  with  6 degrees  of  freedom. 

Here  c^j  denotes  the  element  In  row  1 and  column  J of  the  Inverse  of 
the  covariance  matrix  of  R^,...,Rg  . Levene  and  Wolfowitz  (1944)  present 
expressions  for  this  covariance  matrix  and  for  E(R^)  . For  N » 200,000  one 
has  E(R^)  - 83,320.63,  E(R2)  - 36,664.49,  E(R3)  « 10,556.38,  E(R^)  = 2296.15, 

£(85)  ■ 411.44,  E(Rg)  • 61.79  and  £(87+)  “ 9.09  . Here  R^^  denotes  the  number 

of  runs  of  length  7 or  more.  Although  a statistic  similar  to  (2)  that  Incorporates  < 

87^  can  be  constructed  and  asymptotically  has  a chi-square  distribution  with  7 I 


degrees  of  freedom,  the  small  value  of  E(R7+)  the  present  case 
encouraged  us  to  work  with  (2),  thereby  avoiding  any  discretization  error 
that  Inclusion  of  might  Induce. 

Let  the  superscript  (1)  denote  the  1th  replication  for  a given 
multiplier.  Also  let  = 1 - P(R^^^)  denote  a probability  Integral 

transformation  so  that  has  the  uniform  distribution  on  (0,1).  Then 

‘(O.t]'"*'’*  o.t.l, 

where  I denotes  the  Indicator  function.  Is  an  empirical  cdf.  Figure  1 
shows  Fp(t)  with  n * 100  for  each  multiplier  and  column  1 of  Table  1 
lists  the  Kolmogorov-Smirnov  statistics  0^  = sup|F^(t)  - t|  . Notice 
that  the  test  falls  to  reject  for  each  multiplier  at  the  a = 0.10 
level  but  rejects  multiplier  III  at  the  a = 0.20  level. 

4.  Chi-Square  Test  Statistic 

To  test  H2  we  chose  a chi-square  goodness-of-fit  statistic.  Consider 
K cells  on  the  unit  Interval  each  of  length  1/K  . Let  N|^  denote  the 
number  of  the  n observations  on  a given  replication  that  fall  Into  the 
Interval  ( {k-l)/K,  k/K  ] . Then  for  a specified  K 

C ■ R 'J.1  ‘"k  ■ ■ B tl  "'k  ■ " 

asymptotically  has  a chi-square  distribution  with  K-1  degrees  of  freedom. 

Choosing  K • 2^^  - 4096  Implied  a cell  width  1/K  - 0.000244140625  and  enabled 
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us  to  test  the  first  12  bits  of  . 


Table  1 

Kolmogorov- Smirnov  Test  Results  D 
n = 100  " 


Multiplier 

Runs 

Up  and  Down 
(1) 

■ 

Chi-Square 

(2) 

Serial 

(3) 

Relative 

Execution 

Time 

(4) 

I 

16807 

0.0856 

0.1008 

0.0970 

1 

II 

630360016 

0.0852 

0.0558 

0.0884 

3.03 

III 

1078318381 

0.1068** 

0.0642 

0.0850 

4.42 

IV 

1203248318 

0.0788 

0.0743 

0.0816 

4.82 

y 

397204094 

0.0542 

0.1052 

0.0758 

2.37 

VI 

2027812808 

0.0919 

0.1071** 

0.0673 

7.28 

VII 

1323257245 

0.0903 

0.0926 

0.0723 

5.31 

VIII 

764261123 

»J 

0.0746 

0.1033 

0.1048 

3.46 

a = 

'”•<“100 » “*100)  • 

Source:  Owen  (1962) 

a 

0* 

.10 

.12067 

.20 

.10563 

Figure  2 shows  the  empirical  cdf's  of  C for  the  8 multipliers  and  column  2 
of  Table  1 lists  the  Kolmogorov- Smirnov  statistics.  Again,  the  test  falls  to 
reject  H2  for  each  multiplier  at  the  a ■ 0.10  level;  however.  It  now  rejects 
multiplier  VI  at  a > 0.20  . 
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5.  Serial  Test  Statistic 

Hypothesis  Is  designed  to  detect  nonuni formlty  when  the  are  taken 

In  nonoverlapping  pairs  or  2-tuples.  One  motivation  for  this  testing  arises 
from  the  theoretical  observation  In  Harsaglla  (1968)  that  the  randomness 
of  k-tuples  becomes  more  suspect  as  k Increases.  The  spectral  and  lattice 
tests  support  this  observation.  In  particular,  see  Hoaglln  (1976)  and 
Harsaglla  (1971).  Ideally  one  would  like  to  test  for  the  uniformity  In 
distribution  of  k-tuples  over  the  k-dimenslonal  unit  hypercube.  In  practice, 
such  testing  Is  excessively  expensive,  even  for  k = 2. 

Let  us  divide  the  unit  Interval  Into  K cells,  each  of  width  1/K  . 

Let  Nji^  denote  the  frequency  with  which  falls  Into  the  square 
( ( (j  - 1)/K,  j/K  ],  ( (k  - 1)/K,  k/K  ] ).  For  fixed  K the  quantity 


5 • r . 1 <"ji,  - 


asymptotically  has  a chi-square  distribution  with  - 1.  Suppose  we  had  chosen 
as  before  K ■ 4096.  Then  there  would  be  = 16777216  cells.  To 
guarantee  a mean  of  5 per  cell  under  would  require  N > 80  million 
observations  per  replication  or  over  8 billion  observations  per  multiplier. 

31 

Since  2 - 1 < 4.3  billion  such  a sample  size  Is  not  possible  using  this 

test  procedure.  Because  of  this  demonstrated  excessiveness  we  chose  K = 128 

which  required  16384  cells.  Implied  a cell  width  of  0.0078125  and,  for  N = 200,000, 
2 

a mean  of  n/K  > 12.21  per  cell  under  H^.  This  choice  of  K enabled  us  to 
study  the  first  7 bits  of  the  coordinates  of  . Figure  3 shows  the  empirical 
cdf  of  S for  each  multiplier  and  column  3 of  Table  1 gives  the  Kolmogorov- 
Smlmov  statistics.  Notice  that  no  significance  occurs,  even  for  a - 0.20  . 
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6.  An  Additional  K-S  Type  Statistic 

Although  the  results  In  Table  1 raise  moderate  concern  about  multipliers 
III  and  VI  only,  the  appearance  of  the  empirical  cdf s in  Figures  1,  2 and  3 
raises  broader  concern.  Since  the  empirical  cdf  is  a tied-down  Brownian 
motion  process  [ 4 ],  the  properties  of  such  a process  may  account  for  the 
drift  apparent  in  many  of  the  curves.  As  a further  check  on  the  empirical 
cdf's  we  studied  an  additional  statistic  for  each  empirical  cdf: 


‘n  ■ / ><“  • 


This  quantity  denotes  the  proportion  of  F (t)  that  falls  below  the  45 

n 

degree  line.  Using  results  in  Dwass  (1958}  one  can  show  that  for  given  n 
has  the  uniform  distribution  on  (0,1)  . Presumably,  values  of 
close  to  0 or  1 are  suspect.  Table  2 lists  X^^  in  column  1 for  the  runs 
up  and  down,  chi-square  and  serial  test  statistics  for  the  8 multipliers.  Notice  X^ 
raises  suspicion  about  multipliers  I,  III  and  VI  at  the  a = 0.05  level 
and  about  IV  and  VII  at  the  a « 0.20  level.  In  particular,  the  results 
for  X|^  Indicate  that  the  empirical  cdf's  for  multipliers  I,  III,  IV,  VI 
and  VII  spend  either  more  (or  less)  time  above  (below)  the  45  degree  line 
than  theory  suggests. 


7.  Anderson-Darling  Test 

Although  the  supplementary  statistic  X^  appears  more  discriminating 
than  the  Kolmogorov-Smirnov  statistic,  it  weighs  deviations  equally,  regardless 
of  where  they  occur  in  (0,1).  In  an  effort  to  assign  more  weight  to  deviations 
in  the  tails  of  the  distributions  we  subjected  the  24  empirical  cdf's  to  the 


n 


r 

Table  2 

Additional  Test  Results 


} 


Multiplier 

J 

’‘n 

Y + 
n 

(1) 

(2) 

Runs  Up  and  Down 

I 

0.2433 

0.5208 

11 

.6811 

.6278 

III 

.2924 

1.4979** 

IV 

.7399 

.6693 

V 

.5097 

.41 31 

VI 

.8406 

.8298 

VII 

.0604** 

1.3526 

VIII 

.7251 

1.2748 

Chi-Square 

I 

.2092 

1 .4235** 

II 

.7110 

.3326 

III 

.8204 

.7583 

IV 

.5567 

.3566 

V 

.1341 

1.2414 

VI 

.9828* 

1.5360** 

VII 

.5094 

1.2628 

VIII 

.7165 

.7647 

Serial 

I 

.0207* 

.8981 

II 

.6969 

1.2124 

III 

.9786* 

1 .2242 

IV 

.9372** 

V 

.7272 

.8903 

VI 

.7318 

.2931 

VII 

.2999 

.7091 

VIII 

.2443 

1.2357 

^Significant  for  two  tailed  test  at  a * 0.05  level. 


♦♦Significant  for  two  tailed  test  at  a ■ 0.20  level. 

^Significance  points  were  computed  for  Y using  Anderson  and  Darling's 
expression  for  the  asymptotic  distribution  In  [ 1 ]. 
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Anderson- Darling  test  (1952,  1954).  The  test  statistic  Is 

= n - t]^/t{1  - t)}dt  . 

Since  F„(t)  has  mean  t and  variance  t(l  - t)/n,  Is  the  Integral 

of  sample  mean-square  errors  normalized  by  their  theoretical  mean-square 

errors.  Anderson  and  Darling  (1952)  give  the  asymptotic  distribution  of 

and  Indicate  that  this  limiting  form  Is  approached  rapidly.  Table  2 

lists  the  Y and  raises  suspicion  about  multipliers  I,  III  and  VI  at  the 
n 

o * 0.20  level. 

8.  Execution  Time 

Although  randomness  considerations  principally  determine  a multiplier's 
acceptability,  efficiency  In  execution  also  plays  a role.  This  Is  especially 
true  when  one  regards  several  multipliers  as  equally  good  with  regard  to 
, H2  and  and  must  decide  which  to  use  In  practice.  Column  4 of 
Table  1 lists  relative  execution  times  for  the  8 multipliers,  based  on  runs 
performed  at  the  University  of  North  Carolina  Computation  Center  with  interrupts 
due  to  other  users  eliminated.  The  wide  disparity  In  execution  times  confirms 
a similar  observation  In  Learmonth  (1975).  If  one  plots  the  multiplier  value 
A against  execution  time,  an  approximate  linear  relationship  appears.  This 
may  be  due  to  the  Increasing  number  of  modulo  reductions  that  occur  per 
multiplication  as  A Increases. 

9.  Conclusions  and  Recommendations 

Table  3 presents  summary  test  results.  They  arouse  serious  suspicion 
about  multipliers  III  and  VI  and  suspicion  about  I,  IV  and  VII.  Therefore, 
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a conservative  user  would  select  multiplier  II,  V or  VIII.  Since  II  Is  In 
relatively  common  use,  as  In  the  simulation  programming  language  SIMSCRIPT  II, 
one  may  wish  to  rely  on  this  choice.  However,  Table  1 clearly  shows  that  V 
Is  the  most  efficient  from  the  viewpoint  of  execution  and  statistical  acceptance. 


Table  3 

Summary  Test  Results^ 


Multipl ler 

°n 

^n 

Y 

n 

I 

S* 

C* 

II 

III 

R** 

s* 

R** 

IV 

s** 

V 

VI 

c** 

c* 

C** 

VII 

R** 

VIII 

5 runs  statistic,  C ■ ch1*square  statistic,  S a serial  statistic; 

A single  asterisk  denotes  significance  at  a > 0.05  and  a double  asterisk, 
significance  at  a ■ 0.20 
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This  paper  presents  the  results  of  empirically  testing  8 alternative 
multipliers  for  a multiplicative  congruentlal  generator  with  modulus  23'-!. 
The  LLRANOOM  random  number  package  [11]  uses  one  of  the  multipliers,  the 
simulation  programming  language  SIMSCRIPT  II  uses  a second  and  the  remaining 
six  are  the  best  of  50  candidate  multipliers  studied  by  Hoaglln  (1976)  using 
the  theoretical  spectral  and  lattice  tests.  The  battery  of  tests  fall  to 
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